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Vortices and the entrainment transition in the 2D Kuramoto model 
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We study synchronization in the two-dimensional lattice of coupled phase oscillators with random 
intrinsic frequencies. When the coupling K is larger than a threshold Ke, there is a macroscopic 
cluster of frequency- synchronized osciUators. We explain why the macroscopic cluster disappears at 
Ke- We view the system in terms of vortices, since cluster boundaries are delineated by the motion 
of these topological defects. In the entrained phase {K > Ke), vortices move in fixed paths around 
clusters, while in the unentrained phase [K < Ke), vortices sometimes wander off. These deviant 
vortices are responsible for the disappearance of the macroscopic cluster. The regularity of vortex 
motion is determined by whether clusters behave as single effective oscillators. The unentrained 
phase is also characterized by time-dependent cluster structure and the presence of chaos. Thus, 
the entrainment transition is actually an order-chaos transition. We present an analytical argument 
for the scaling Ke ~ Kl for small lattices, where Kl is the threshold for phase-locking. By also 
deriving the scaling Kl logN, we thus show that Ke ^ logN for small A'^, in agreement with 
numerics. In addition, we show how to use the hnearized model to predict where vortices are 
generated. 

PACS numbers: 



I. INTRODUCTION 

Collective behavior of coupled oscillators is found in 
many areas of science [E Examples include Joseph- 
son junctions 0, lasers [3], neural networks chemical 
oscillators and nanomechanical resonators Re- 
cently, there has been much interest in studying coupled 
oscillators on complex networks, motivated by biological 
and social networks [sl-fTTj. 

This field is an interesting marriage of statisti- 
cal physics and nonlinear science, because these non- 
equilibrium systems may exhibit phase transitions. A 
common approach is to consider populations of oscilla- 
tors with random intrinsic frequencies. The coupling 
between oscillators acts against the frequency disorder 
to synchronize the oscillators. The original Kuramoto 
model considered oscillators with all-to-all coupling Q. 
The model has also been studied on low-dimensional lat- 
tices with local interactions [T2| - [l6| and long-range inter- 
actions [IMl- 

Two synchronization transitions occur on low- 
dimensional models as the coupling K changes. The en- 
trainment transition at Ke denotes the onset of macro- 
scopic entrainment, when there is a cluster of frequency- 
synchronized oscillators on the order of the system size. 
There is also a transition at Kl to the phase-locked state, 
in which all oscillators evolve with the same frequency. 

In the tradition of statistical physics, the main question 
is whether or not the entrainment transition exists, i.e., 
whether Ke is finite, in the limit of infinite system size 
{N — >■ oo). Simulations indicate that in the case of local 
interactions, it exists only in dimension d > 3, meaning 
that the lower critical dimension is 2 [l^. Ideally, one 
would show this analytically, but it is difficult due to 
the presence of nonlinearity and disorder. It has been 



proven that macroscopic entrainment does not exist in 
d = 1 in an infinite system [ij, [ij]. Various heuristic 
arguments indicate that it exists only in d > 3 p^[T3l . [20| . 
However, there has been no clear explanation of exactly 
what happens at the entrainment transition. 

In this paper, we examine the entrainment transition 
in the 2D model and elucidate why and how it happens 
in a finite system with local interactions. Surprisingly, 
there are several differences between the entrained phase 
{K > Ke) and the unentrained phase (K < Ke) besides 
the presence of the macroscopic cluster. The entrained 
phase is characterized by time-independent cluster struc- 
ture, while cluster boundaries continually change in the 
unentrained phase. In fact, the entrained phase is not 
chaotic while the unentrained phase is chaotic, meaning 
that the entrainment transition is actually an order-chaos 
transition. 

It is convenient to view the system in terms of vortices, 
since the boundaries of frequency clusters are delineated 
by vortex paths. In the entrained phase, vortices move 
in fixed paths around clusters, while in the unentrained 
phase, vortices sometimes wander off. These deviant vor- 
tices are responsible for the system- wide detrainment. 

In an A^ X A^ lattice, the nature of the entrainment tran- 
sition is different in small (A^ < 50) and large (A^ > 100) 
lattices. In small lattices, the transition is determined by 
a cluster made by a vortex pair, while in large lattices, the 
transition is due to clusters made by single vortices. We 
focus on small lattices and provide analytical arguments 
for the scalings Ke ~ Kl and Kl ^ logA^. Hence, we 
derive the scaling Ke ~ log A^ for small A^, in agreement 
with numerics. We also show how to use the linearized 
model to predict where vortices are created. 

Previous works have touched on certain aspects of our 
results. Topological defects have been observed in 2D 
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lattices of chaotic oscillators (where individual oscilla- 
tors are intrinsically chaotic) [2l[ . Chaos and phase slips 
have been observed in the ID Kuramoto model [131 ■ Fre- 
quency clusters have been studied also in the ID Ku- 
ramoto model [2^ [13] ■ This paper explains how all these 
ideas are related to each other and to the entrainment 
transition. 

The outline of the paper is as follows. In Sec. II, we 
review the 2D model. In Sec. Ill, we summarize qualita- 
tively its complex phenomenology. Then we present an- 
alytical results on cluster stability in Sec. IV. We study 
phase-locking in Sec. V and conclude in Sec. VI. 

II. MODEL 



Although this quantity diverges as — 00, it does not 
rule out the possibility of entrainment in 2D; there may 
be frequency order in the presence of large phase gradi- 
ents, when the linear model is not applicable. Also, one 
might guess incorrectly from Eq. ^ that Ke and Kl 
scale as (logA^)^. 

The average frequency of an oscillator is defined as 
^ = [(^{to+T) — 0{tQ)]/T, where to and T are the transient 
and averaging times, respectively. An oscillator is said 
to be frequency-synchronized with its neighbor if their 
phase difference shifts by less than tt during the averaging 
time. In this paper, the numerical integration of Eqs. ([T]) 
was done using the Euler method with to = T = 10*, a 
time step of 0.02, and initial phases set to zero. 



We consider the N x N two-dimensional lattice of os- 
cillators with nearest-neighbor interactions and periodic 
boundary conditions: 
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The intrinsic frequencies are Gaussian distributed 
with zero mean and unit variance. We assume, with- 
out loss of generality, that the average of w in a given 
realization is zero. 

We also consider the linearized version: 

Oij = uJij + K[{9i^ij — 9ij) + — 9ij) 

+ %) + (%+!-%)] (2) 

This approximation is sometimes useful, because for large 
enough K and finite A^, the lattice has small phase 
gradients and is thus in the linear regime. The ob- 
vious advantage of the linearized model is that it is 
straightforward to solve via discrete Fourier transforms. 
With the definitions Oki = Emn and 
Wm„e"'T^(''™+'"\ each 4; satisfies 
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Each Fourier component decays exponentially towards its 
steady-state value so that at steady state, 9mn = (^mln 
where 
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From Eq. ©, {6*^^';} also satisfy: 
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The standard deviation of differences between neighbor- 
ing phases can be calculated in the continuum approxi 
mation as UM- 



1 /logA^ 
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III. ENTRAINED AND UNENTRAINED 
PHASES 

The 2D Kuramoto model has a rich phenomenology. 
In this section, we qualitatively describe the behavior in 
order to motivate analytical calculations in later sections. 



A. Clusters 

Consider what happens as K changes. Above ATl, all 
oscillators have the same average frequency. Immediately 
below Kl, & small cluster of frequency-synchronized os- 
cillators appears, while the rest of the system constitutes 
a macroscopic cluster [Fig.[IJa)]. As K decreases further, 
more small clusters appear, and they generally maintain 
the same shape as K decreases. Below Ke, there is no 
longer a macroscopic cluster. 

There are three qualitative differences between the en- 
trained and unentrained phases besides the presence of 
a macroscopic cluster. The first difference is the con- 
stancy of clusters. Above Ke, the cluster structure is 
time-independent after a sufficient transient time. How- 
ever, below Ke, the cluster structure changes over time, 
so it is hard to say which oscillator is synchronized with 
which, since it depends on the values of to and T. This 
time-dependence is surprising, since one would expect the 
unentrained phase to have small but well-defined clusters. 
(For example, the ID chain with random uj and K has 
a time-independent cluster structure, even in the unen- 
trained phase [23l. \2^). 

The second difference is in the distribution of average 
frequency differences between neighbors (Fig. [5]). Differ- 
ences in average frequency are due to 2n phase slips. In 
the entrained phase, a pair of neighbors will experience 
either zero or many phase slips during the averaging time. 
But in the unentrained phase, some pairs experience only 
one or few phase slips. For a given disorder realization, 
as K is decreased below Ke, there is a sudden drop in 
the size of the largest cluster and a sudden appearance 
of single 27r phase slips (Fig. [3]). 
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FIG. 1: Example of a 50x50 lattice in the entrained phase, (a) 
Average frequency, showing the two frequency clusters, (b) 
Time-averaged vorticity, where black and white are opposite 
polarities. Vortex paths delineate the cluster boundaries, (c) 
Phase (mod 2ti) at a point in time, where black is and white 
is 2-n. (d) The equivalent lattice of vortices at that time, where 
black and white are opposite polarities. In (a,c), each pixel 
represents an oscillator. In (b,d), each pixel represents a unit 
cell of four oscillators. 
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FIG. 2: Distribution of the number of 2% phase slips between 
neighboring oscillators for a 50 x 50 lattice, (a) In the en- 
trained phase, neighbors have either zero or many phase slips 
during the averaging time, (b) In the unentrained phase, 
some neighbors have single or few phase slips, (c) and (d) 
are zoomed-in views of (a) and (b), respectively. Here, the 
entrainment transition happens at Ke = 1.32. 




The third difference is that the entrained phase is not 
chaotic while the unentrained phase is. The largest Lya- 
punov exponent is zero in the entrained phase but greater 
than zero in the unentrained phase . For a given lat- 
tice, the onset of chaos occurs at the same K as the 
appearance of single phase slips (Fig. The entrain- 
ment transition usually happens at the same K , although 
sometimes lower. 

Averaging over the disorder. Fig. 2] shows that Ke ^ 
logN, where Ke is defined as when the largest cluster 
encorrrpasses half the lattice. This scaling agrees with 
Ref . [l^ . The figure also shows that the onset of chaos 
and single phase slips coincides well with the entrainment 
transition. 
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FIG. 3: Size of the largest cluster, number of single phase 
slips, and largest Lyapunov exponent vs. coupling K for a 
50 X 50 lattice. In this disorder realization, the entrainment 
transition happens at Ke = 1.48, and phase locking happens 
at Kl = 1.58. As K decreases past Ke, the contraction of the 
largest cluster coincides with the appearance of single phase 
slips and chaos. 



These observations indicate that the system-wide de- 
trainment is caused by the random propagation of single 
phase slips. Indeed, Fig.[5]shows that when K < Ke, the 
largest cluster shrinks over time due to occasional phase 
slips within it, cutting it up until there is no longer a 
macroscopic cluster. 



B. Vortices 

It is useful to speak of vortices to describe the propa- 
gation of phase slips. Looking at the numerical evolution 
of the solution, one sees that the phase is mostly smooth 
except for small regions where the phase gradient is large 
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FIG. 4: Coupling threshold for phase-locking Kl (black tri- 
angles), entrainment Ke (blue circles), and onset of chaos 
and single phase slips (red squares) for lattices with A'^^ os- 
cillators. Each data point is averaged over 50 disorder real- 
izations. The standard deviation of the mean is about 0.02 
for all points. Both Kl and Ke scale as log A^. The dashed 
line plots 0.44 log A'^. Here, Ke is defined as when a cluster 
encompasses half the lattice. 
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where the phase differences have been shifted mod 27r 
into the range (— tt, tt]. The curl can be -|-27r, — 27r, or 
0, corresponding to a -I- vortex, a — vortex, or no vortex 
[Fig. HTd)]. Note that a vortex exists on a unit cell of 
four oscillators, not just on one oscillator. 

From the above definition of a vortex on a discrete 
lattice, it follows that there must be an equal number 
of -|- and — vortices in the case of periodic boundary 
conditions. Also, the curl on a unit cell changes only 
when one of its crosses tt (mod 27r). Thus, a vortex 
moves to a neighboring cell depending on which A0 slips. 
In other words, vortex motion is equivalent to phase-slip 
propagation. Lastly, vortices are created when a pair of 
oscillators without neighboring vortices slips, resulting 
in two vortices of opposite polarity, one on either side. 
There are certain places in a lattice that tend to create 
vortices. 

Since a vortex is topological, it exists until it meets 
and annihilates with a vortex of opposite polarity. Vor- 
tex paths delineate the boundaries of frequency clusters 
[Figs. [TJa) and[Hb)]. When a vortex passes between a 
pair of oscillators, it causes a 27r phase slip. The accu- 
mulation of phase slips, due to repeated vortex crossings, 
leads to differences in average frequency. 

Thus the entrained phase is characterized by vortices 
moving in fixed paths around the clusters, meaning that 
the vortices are locally confined and the cluster bound- 
aries are time- independent. In the unentrained phase, 
vortices move inconsistently and sometimes wander off 
deeply into the formerly macroscopic cluster and chop it 
up (Fig. [5]). A signature of inconsistent motion is the 
presence of single phase slips during the averaging time, 
due to vortices that passed by only once, as opposed to 
regularly [Fig.[2Ud)]. This irregularity reflects the chaotic 
nature of the unentrained phase. 

In the rest of the paper, we call a cluster stable when its 
vortices move consistently and unstable when they move 
inconsistently. Instability connotes the presence of single 
phase slips and chaos. The fact that the entrainment 
transition occurs at a slightly lower K than the onset 
of instability (Fig. indicates that the lattice becomes 
unentrained because the clusters become unstable. 



FIG. 5: Relative size of the largest cluster vs. averaging time 
T for three different 50 x 50 disorder realizations. K is below 
Ke, so the lattices are in the unentrained phase. The down- 
ward jumps are due to occasional wayward vortices that split 
up the cluster. This also shows that it is necessary for T to 
be large in order to accurately distinguish the entrained and 
unentrained phases. 



[Fig- [He)]. These are vortices, which are topological de- 
fects arising from the 27r-periodicity of each 9. The phase 
winds by 27r around a vortex. To identify vortices, one 
computes the lattice curl of the phase gradient around 
each unit cell of four oscillators [26] . The curl is equal to 
the sum of the directed phase differences around a cell. 



C. Two types of entrainment transitions 

As K decreases below Ke, the lattice transitions from 
the entrained to the unentrained phase. Simulations indi- 
cate that small lattices {N < 50) have a different pathway 
to the unentrained phase than large lattices {N > 100). 
The difference between the two is in the nature of the first 
microscopic cluster that appears at the locking threshold 
Kl. In small lattices, that cluster is usually made by 
a vortex pair, while in large lattices, it is usually made 
by a single vortex. Below, we describe the two types of 
clusters and their corresponding entrainment transitions. 
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1. Vortex-pair cluster 

The first type of cluster is made by a pair of oppositely- 
charged vortices that are created at a certain spot, travel 
along the cluster boundary, and then annihilate with 
each other. This cycle repeats periodically over time, 
so that the cluster has a different average frequency from 
its neighbor. We call this a vortex-pair cluster^ and an 
example is in Fig. [T] 

As K decreases, the shape of the cluster stays the 
same, because the vortices travel along the same paths 
although the vortices are produced more and more fre- 
quently. When K < Ke, the cluster is unstable and pro- 
duces vortices which sometimes deviate from the original 
path. The liberated vortices sometimes cut across the 
macroscopic cluster so that the system becomes unen- 
traincd. Since a vortex exists until it annihilates with 
an oppositely charged one, a deviant vortex may travel a 
long distance before annihilating. 

The stability of a vortex-pair cluster (and thus the 
existence of the entrained phase) is determined by the 
balance of two time scales: the lifetime of a vortex pair 
(tiife) and the period of pair production (tper )- When the 
cluster is formed at K^, tufe <C tper- As K decreases, 
tiife increases while tper decreases (Fig. [6]). The cluster 
becomes unstable when tuje ~ tper- In fact, a cluster 
with tiife > tper IS ncvcr stable. This means that the 
entrainment transition happens when a new pair is pro- 
duced immediately after the previous pair annihilates. 
So it seems that the transition is due to the balance of 
time scales, instead of the interaction between vortices 
of different clusters. We discuss this case further in Sec. 

m 



2. Single-vortex cluster 

Another type of cluster is made by a single vortex that 
continually orbits it. We call this a single-vortex cluster. 
It is usually composed of only one or a few oscillators. 
In the case of periodic boundary conditions, topologi- 
cal constraints require that there be an even number of 
vortices, so there are an even number of single-vortex 
clusters (Fig.[7|). 

The stability of single-vortex clusters is an interest- 
ing phenomenon that invokes the topological nature of a 
vortex. Suppose K > Kl, so that all the oscillators are 
phase-locked. Due to the disorder of intrinsic frequen- 
cies, there is a time-independent phase gradient across 
the lattice. As K decreases, the phase differences be- 
tween neighbors increase. When K ~ Kl, there is a pair 
of oscillators on the verge of slipping relative to each 
other, and their phase difference is the largest in the lat- 
tice. When K decreases below Kl, that pair does slip, 
producing two vortices that move. But instead of meet- 
ing and annihilating, the vortices are each pinned to a 
single-vortex cluster. A vortex has a phase field around 
it that accumulates by 27r. The oscillator pair that pro- 
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FIG. 6: Ratio of vortex lifetime to the period of vortex pro- 
duction for the cluster formed at Kl ■ Each line corresponds to 
a disorder realization of size N xN and with a cluster contain- 
ing m oscillators. From left to right, (N,m) is (10,3), (10,1), 
(10,1), (25,64), (25,60), (50,14), (25,6), and (50,260). The ra- 
tio is when the cluster is formed at Kl and increases as K 
decreases until the ratio is on the order of 1. At that point, the 
cluster is unstable, and the vortices are sometimes liberated 
from their orbit around the cluster, leading to system-wide 
detrainment. 

duced them is actually prevented from slipping again by 
the phase field of the vortices: a vortex is positioned rela- 
tive to the oscillator pair, so that its accumulating phase 
field balances the tendency of the pair to slip. Thus no 
more vortices are produced, and the cluster configuration 
is stable. 

As K decreases further with the vortices in this con- 
figuration, the phase difference of the oscillator pair in- 
creases again. When K is low enough, it finally slips and 
produces more vortices. If the vortices do not find a new 
configuration to stop the oscillator pair, it will continue 
to produce vortices that detrain the lattice. Therefore, 
whether the lattice is entrained depends on whether the 
source of vortices is quenched. When it is quenched, the 
vortices are locally confined to single- vortex clusters, but 
when it is not quenched, the vortices wander throughout 
the system. 

IV. STABILITY OF VORTEX-PAIR CLUSTER 

Having described the features of the 2D model, we now 
do some analytical calculations. In particular, we are 
interested in how the entrainment threshold Ke scales 
with the system size N . Since, a lattice becomes unen- 
trained when the microscopic cluster formed at the lock- 
ing threshold Kl becomes unstable, we study the onset of 
instability. Here, we calculate Ke{N) for small N, where 
the entrainment transition is determined by a vortex-pair 
cluster. We leave the stability of single-vortex clusters 
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FIG. 7: Example of a 100 x 100 lattice with two single-vortex 
clusters. Only a portion of the lattice is seen here, (a) Av- 
erage frequency, showing the two single-vortex clusters and 
the macroscopic cluster, (b) Time- averaged vorticity, where 
black and white are opposite polarities. The vortices are each 
pinned to a cluster, (c) Phase (mod 2n) at a point in time, 
where black is and white is 27r. (d) The equivalent lattice of 
vortices at that time, where black and white are opposite po- 
larities. The arrow points to the source of the vortices, which 
does not produce any more vortices. 



(and hence Ke{N) for large N) for future work. 

Our approach is based on the fact that the unentrained 
phase is chaotic while the entrained phase is not. Con- 
sider the interaction of the microscopic cluster with the 
macroscopic cluster that surrounds it. When the two 
clusters are stable, each may be considered as a single ef- 
fective oscillator, since the constituent oscillators within 
each behave coherently. A system of two coupled os- 
cillators is not chaotic, which is consistent with the as- 
sumption that the clusters are stable. However, when 
the clusters stop behaving as effective oscillators, the new 
degrees of freedom within each cause the system to be- 
come chaotic. In the presence of chaos, the vortices are 
no longer confined and proliferate to cause system-wide 
detrainment. 

So the question is: when does a cluster stop behaving 
as an effective oscillator? Recall that a vortex-pair cluster 
is made by periodically produced vortex pairs. There are 
two time scales involved: tu fe is the lifetime of the vortex 
pair and iper is the period of pair production. In other 
words, tiife is how long it takes for a pair to annihilate, 
and iper is the time between the creation of successive 
pairs, tiiff, is also the duration of the cluster- wide phase 
slip. By definition, an effective oscillator can have only 



one phase slip during iper- Thus, a ckister cannot behave 
as an effective oscillator when tufe > tper- Thus the 
stability condition is 

When this is satisfied, the clusters behave as effective 
oscillators [l^. We expect this condition to be satisfied 
for large but not small K. Note that this explanation is 
consistent with simulation results (Fig. [B]). 

We proceed with a self-consistent argument: assum- 
ing that Eq. ([7]) is satisfied, we calculate tufe and tp^r 
and then see when the condition is no longer satisfied. 
We consider the interaction between the first microscopic 
cluster and the macroscopic cluster that surrounds it. 
Simulations indicate that as K decreases, the shape of 
the first cluster generally remains unchanged. Hence, we 
assume here that the first cluster keeps the same shape 
for K > Ke- 

Let the first microscopic cluster and the macroscopic 
cluster be denoted by A and B, respectively. We coarse- 
grain the niA oscillators of cluster A into a single effec- 
tive oscillator with phase 9^ — X^^gA intrinsic 
frequency oj^ — J2ieA^i- Let 69^ describe an oscil- 
lator's deviation from the cluster phase: 9i = 9"^ + 69 f. 
We define similar quantities for cluster B. Since Eq. ([7]) 
is assumed to be satisfied, S9f and 69f can be taken 
to be time-independent. In other words, the phases are 
rigid within each cluster. 

The phase difference tf) = 9^ — 9^ satisfies 

= AOJ--Y, M<P + S9f - 69 f) , (8) 

where Aw = wa — and u = "^^"^^ . The sum runs 
over the lattice edges that connect A and B, since the 
coupling terms within each cluster cancel due to action- 
reaction symmetry. Let q be the number of edges between 
A and B. Then we write the sum of sines as a single sine: 

aK 

(j> = Aoj sin(0 + /3) , (9) 

where 

a{K)' = q+ Yl cos{{69f-69f)-{69t-69f)), 

(10) 

and (3{K) is a phase offset that does not matter in the fol- 
lowing discussion. The deviations 69i depend on K and 
can be estimated using the linear solution (Sec. IVB II) . 
In the limit of large K, all 69i sa 0: the coupling terms 
add coherently and a — q. When K is small, 69 i is large: 
the sum in Eq. ()10p is over q{q — 1) random numbers so 
y/q. 

The period of vortex production tp^r is the period of 
(j), since each vortex pair passing between A and B con- 
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(11) 

(12) 



where r — 



When tper diverges, cluster A is syn- 



M|Aw| • '■'^i'^ii '•per 

chronized with cluster _B and the lattice is phase-locked. 
Thus we identify ulK^ = hAuj. Assuming that a ~ aL 
in the range of K that we are interested in, 



iper — 



27r/i 



(13) 
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X as in a 



Note that when K Kl, 
one-dimensional ring . 

We now estimate the lifetime of a vortex pair tu fe ■ The 
presence of a vortex pair is a manifestation of the fact 
that the clusters are experiencing a phase slip relative to 
each other. Hence, we should calculate how long it takes 
the clusters to slip. First, we calculate the duration of 
a phase slip tps in the effective-oscillator model. This is 
not the same as i/i/e, since vortices are not present in 
the effective-oscillator model. We will later add in the 
vortices by accounting for the fact that each cluster is 
spatially distributed. To estimate tps, we calculate the 
duration of time when is large: 



1 



d<j) 
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TT — 2 tan 
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(14) 



(15) 



For the sake of scaling, the exact limits on the integral 
do not matter. As r decreases from 1 to 0, tps increases 
slightly from ^ to Thus, 



1 A* 



\Auj\ ulKl 



(16) 



To convert tps to tu fe , we use the fact that a vortex pair 
is created when the first edge phase difference between 
A and B crosses tt (mod 27r) and annihilates when the 
last edge phase difference crosses tt (Sec. IIIIBI) . Thus we 
consider 



7(if) = maMSei 



9f) - mHsei 



9f), (17) 



where each term corresponds to one of the edges between 
A and B. When a vortex pair is produced, 7 is how far 
(j) needs to go before the pair annihilates. Due to the dis- 
order of intrinsic frequencies, 7 increases as K decreases. 
Since (0) ^ j— during the lifetime of a vortex pair, 



Fortunately, it is not necessary to calculate ^(K) ex- 
plicitly. Since is defined as when tuf^. « ipen we 
know that ^(Ke) ~ 27r, because then a vortex pair is 
created immediately after the previous one annihilates. 
We find Ke hy equating Eqs. (IT3)) and (fT8|) and arrive at 
the scaling 



Kr 



(19) 



Thus, the entrainment transition is tied to the phase- 
locking transition. In Sec. |Vl we show that Kl ~ log N . 
This means that Ke logiV, in good agreement with 

Fig.m 

Note that our argument only applies to lattices where 
the first microscopic cluster is a vortex-pair cluster [N < 
50). We have not analytically studied the stability of 
single- vortex clusters, but Fig. |4] indicates that Eq. (fT9| 
would also apply to those cases. 



V. PHASE-LOCKING 

For a given lattice of oscillators, when K > Kl, the 
system is phase-locked and all the oscillators have the 
same frequency. In this section, we calculate how K^ de- 
pends on system size N. Since K^ is different for different 
realizations of the intrinsic frequencies, we are actually 
interested in how the disorder-averaged (Ke) depends on 
N. We first consider ID and then 2D. Although phase- 
locking in ID has already been solved we review 
it and then redo it using a linear approach in order to 
tackle 2D, where the usual approach does not work. 



ID 



1. Usual approach 

Consider a one-dimensional chain of N oscillators with 
open boundary conditions: 

61 = uji + Ks'm{02 - 9i) 

e, = w, + X[sin(6ii_i-6l,)+sin(6ii+i-6'i)] 

i = 2,...,N ~l (20) 
9n = lom + K su\{9n-i — 9m) 

We assume, without loss of generality, that the average lo 
is 0. Then the phase- locked solution is given by 9i, such 
that all 9i = 0. Due to the open boundary conditions 
and action-reaction symmetry, one can solve for A^^ = 
9i+i — 9i by adding up the first i equations: 



tlife ~ l{K) tps 



ulKl 



(18) 



K sin A9, 



l<i<N 



(21) 
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A necessary and sufBcient condition for the existence of 
a phase-locked sohition is that: 

i 

max{|^a;,|}<if. (22) 

J = l 

Thus, Kl = max{| "^il}- There are 2^^^ sohitions 

of Eq. (mi, since each A0i can have two values. 

The unique stable solution is given by all cos A6'i > 0. 
One can show this by considering the Jacobian Jij = ^ 
at the phase- locked solution and requiring ^j^- XiJijXj < 
for all perturbations {xi}, such that ||x|| > and 
Y^7=i — ^- The last condition is due to the fact 
that there is always a eigenvalue of J with eigenvector 
(1, 1, • • • , 1)"^, corresponding to the uniform displacement 
of all phases. One finds: 

^^XiJijXj -K ^ {xi - Xi+i)"^ cos A9i . (23) 

ij l<i<N 

Due to the arbitrariness of {xi}, a necessary and sufficient 
condition for stability is that all cos A9i > 0. So only one 
of the 2^~^ phase-locked solutions is stable. 

Using Eq. P^ . it is possible to derive the scaling 
Kl ^ Vn in ID based on random- walk arguments [T5| . 



Eqs. (p4|) with the original nonlinear term. In order to 
get the same amount of coupling force to maintain ^ = 0, 
the corresponding A.9 must be increased in magnitude. 
This is seen in Eq. ([25]). However, if |A6''™| > 1, it is 
impossible to get the same amount of coupling force, so 
there is no equivalent phase-locked solution. 

Consider the dynamics. Suppose the nonlinear model 
starts with all 9i = 0. Since the phase differences are 
small, the dynamics are initially linear, and the phases 
approach the linear steady state. If |A0-*"| < 1, A9i 
converges to sin'^ A6'^*", but if |A6'f"| > 1, A6'i cannot 
converge so it continues to increase and phase slips. Af- 
ter the phase slip, all phase differences are small again 
(mod 27r), so the dynamics are linear for sometime until 
A6', slips again. Thus, if is large so that all |A0""| < 1, 
the system is phase- locked. If one |A0'™| > 1, there 
will be a break in the chain there. If more than one 
|A6''*"| > 1 and their density is small, each one corre- 
sponds to a break. 

Thus, we are doing a local self-consistency check on the 
linear solution. If the linear solution has a big phase gra- 
dient somewhere, it is locally inconsistent there and will 
phase slip. We also note that this is an easy way to find 
the phase-locked solution in ID with periodic boundary 
conditions. 



2. Linear approach 

Consider the linearized version of Eqs. ([20)1 : 

9i ^ LJi + K{92~9i) 

9, = uj, + K[{9,^i - 9i) + {9,+i - 9,)] 

i = 2,...,7V-l (24) 
On — + K{Sn-i — On) 

The solution to the linear model is straightfoward to find. 
Let the steady state of the linear model be {^f "}, which 
always exists. One may get the phase-locked solution of 
the corresponding nonlinear model by making the ansatz: 

AO^ = siij-^ A6if" . (25) 

Plugging this into Eqs. ([20| . one immediately returns to 
Eqs. (1^^ . which are all zero in the steady state. Thus, 
it is easy to go from the linear solution to the nonlinear 
phase-locked state. In light of this, a necessary and suffi- 
cient condition for phase-locking in the nonlinear model 
is 

max{|A6'f"|} < 1 , (26) 

which is equivalent to Eq. (P^. 

This approach has a nice intuitive interpretation. The 
coupling force is stronger in the linear model, since A9 
is steeper than sin A9. Imagine starting from the lin- 
ear steady state and replacing a single linear term in 



B. 2D 

1. Phase-locked solution 

In 2D, the usual approach does not work for finding 
the phase-locked solution, either with open or periodic 
boundary conditions. This is because the connectivity 
prevents one from isolating a single sine term like in 
Eq. ((2T|) . However, the linear approach is still applicable 
in 2D, since the linear solution, Eq. always exists. 
In the next section, we use this approach to show that 
KL-\ogN. 

There is a catch though in 2D due to the presence of 
loops. Suppose one has the linear steady-state solution 
{6*^}"}. The linear approach says that A9 = sin"^ A9"'' 
for each horizontal and vertical pair of oscillators. This is 
actually not a physically allowed solution since the sum of 
A9 around a unit cell would generally not be a multiple of 
27r. This was not a problem in ID because there were no 
loops. Thus in 2D, the linear approach is no longer exact. 
However, Fig. [5] shows that the relationship between the 
actual A9 and A^'*" is still approximately an inverse sine. 

There is another difference in 2D due to loops: it is pos- 
sible to have a stable phase-locked solution with cos A9 < 

0. Consider the Jacobian of Eqs. ([!]): Jijki — The 
stability condition is that J^iju^ij'^ijki^ki < for all 
perturbations {xij}, such that ||x|| > and ^ - ■ Xij — 0. 
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FIG. 8: The actual phase-locked solution vs. the linear steady-state solution for (a) K = 1.2Kl and (b) K — Kl- 

There is an approximate inverse-sine relationship (cyan line). In (b), the phase-locked solution sometimes has \A0\ > ^ due to 
the presence of loops in 2D. Both plots are based on the same 50 x 50 lattice. 



One finds 

xJx = -K^^[{xij - Xi-ij)'^ cos{6i^ij - 9ij) 

+ {xi.j - Xij^if cos(%_i - 9ij)] . (27) 

One might think that a necessary condition for stability 
is that all cos > 0, since if a cos(0i_ij — %) < and all 
x-differences but Xij — Xi-ij were 0, then xJx > 0. The 
catch is that one cannot set all but one x-difFerence to 0; 
due to the loop nature there must be at least one other 
nonzero a;-difFerence to compensate. Thus the require- 
ment that all cos > is only a sufficient condition for 
stability. Indeed, for K « Kl, the stable phase-locked 
solution has some cos/S.9 < [Fig. IHJb)]. 

Note that these are features not just of the 2D lat- 
tice, but of any oscillator network with loops. In 2D, 
the core issue is that there are phases but 2iV^ phase 
differences, so specifying the differences overspecifies the 
system. 

2. Phase-locking cntenon 

In ID, there was phase locking if and only if 
max |A6''™| < 1. This raises the question whether there 
is also a critical value A6'^™ for max |A6'""| in 2D. The 
critical value is probably not an absolute constant since 
the linear approach is an approximation in 2D. Accord- 
ing to simulations, the critical value is actually narrowly 
distributed around 1.27 and seems to be independent of 
system size (Table |l|. The fact that the value is fairly 
consistent means that the intuition behind the linear ap- 
proach is still valid: if the linear solution has small phase 
gradients everywhere, then there is phase-locking, but 



TABLE I: Critical values of phase differences of the linear 
steady-state solution ior N x N lattices. If max |A6''™| < 
A9'^" , the system is phase-locked; otherwise, there are vortex 
sources. These values were found numerically for 50 realiza- 
tions for each size. A9'^" is randomly distributed with the 
sample standard deviation given by the digit in parantheses. 



TV 




10 


1.28(7) 


25 


1.28(6) 


50 


1.26(5) 


100 


1.27(5) 


200 


1.25(5) 



if it has a big gradient somewhere, that place will be 
a source of vortices. Thus, the linear solution predicts 
where the first vortices are formed. 

This provides a way to calculate the critical coupling 
for phase- locking Kl, since A6''*" is inversely related to 
K in Eq. (j4|). This reasoning may seem circular, since 
we first found A0^™ by numerically checking for phase- 
locking and are now using it to determine the phase- 
locking criterion. Actually, one can derive the scaling 
Kl{N) just by positing the existence of some A0^™. The 
actual value for A6'J,™ is only in the prefactor, as seen in 

Eq. dsni). 

We first calculate the probability of phase-locking 
P{K, N) for a given coupling K and system size N. The 
condition for phase- locking is that all |A6''*"| < A6'^*". 
Since is a random variable with standard devia- 

tion a^giin given in Eq. (jH]), P{K,N) is the probability 
that all 2iV^ phase differences are less than the threshold. 
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Assuming that the A6 are independent 



VI. CONCLUSION 



P{K,N) 



erf I Ml"'K 



27r 
loglV 



-, 2N^ 



(28) 



The distribution of Kl is Pk^K) = ^^W^- The 
disorder-averaged Kl is estimated by the K that maxi- 
mizes ■ One gets the impUcit equation 



erf (^A6'^*"(ifi 



2tt 
logN 



{Kl 



2N^ 



1 /log TV 



8tt 



■ exp 



27r(A6'f')2(i^i)2 
iog]v 



(29) 



Taking the logarithm of both sides and considering the 
largest terms in the limit of large TV, 



{Kl 



1 



AogN K 0.44 log TV 



(30) 



which agrees well with numerics (Fig.|4]). Thus, Kl scales 
as log TV. 

The scaling JC r, ~ log TV is consistent with bounds in 
a previous work [l^. A necessary condition for phase- 
locking is that all jwij < AK . The probability of phase- 
locking is bigger than the probability that all TV^ w's 
satisfy this, leading to a lower bound on Kl that scales 
as VlogTV. 

We also note that the above approach to calculate Kl 
is appropriate when the initial phases are zero, since the 
dynamics are linear. When the initial phases are ran- 
domized, the system might not lock when K > Kl- The 
clusters in this case are due to vortices created by the 
initial conditions instead of spontaneously. 



In conclusion, we have studied the entrainment and 
phase-locking transitions in the 2D Kuramoto model. 
We derived the scaling Ke ~ log TV for small lattices, 
in agreement with simulations. We have relied on the in- 
sight that in 2D, a system is unentrained when a cluster 
stops behaving as an effective oscillator, allowing vortices 
to be liberated. 

The next step would be to analytically derive Ke{N) 
for large lattices, since one would like to know whether it 
is finite in the limit TV — >■ cx). This would require explain- 
ing the stability of single-vortex clusters, which seems 
like a challenging theoretical problem. However, doing so 
would answer an open question regarding the Kuramoto 
model. It would also be an interesting problem to pre- 
dict cluster structure using a real-space renormalization 
group approach, as was done in ID [23, 24]. 

Finally, one should study the role of topological de- 
fects in other synchronization models, such as higher di- 
mensions, complex networks, long-range interactions, or 
other oscillator types. For example, one can create a 
small- world network by randomly rewiring some links in 
a 2D lattice, so that the average distance between oscilla- 
tors drops dramatically [11] . One should see whether the 
unentrained phase is still characterized by inconsistent 
vortex motion, chaos, and time-dependent cluster struc- 
ture. Thus, the results in this work may prove useful in 
building a solid mathematical understanding of synchro- 
nization in a variety of situations. 

This work was supported by Boeing. GR thanks the 
Research Corporation and the Packard Foundation for 
their generous support. 
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